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1. Motivation 

In a previous report (Jimenez et al. 1993, referred from now on as JWSR), 
it was shown that the strong vorticity in isotropic turbulence is organized into 
tubular vortices (“worms”) whose properties were characterized through the use 
of full numerical simulations at several Reynolds numbers. That paper should be 
consulted for details as well as for the previous history of the subject. At the time 
most of the observations were kinematic, and several scaling laws were discovered 
for which there was no theoretical explanation. In the meantime, further analysis 
of the same fields yielded new information on the generation of the vortices, and it 
was realized that even if they have to had to be formed by stretching, they were at 
any given moment actually compressed at many points of their axes (Jimenez and 
Wray, 1994). This apparent contradiction was partially explained by postulating 
axial inertial waves induced by the nonuniformity of the vortex cores, which helped 
to “spread” the axial strain and allowed the vortices to remain compact even if not 
uniformly stretched. The existence of such solutions was recently proved numerically 
by (Verzicco, Jimenez & Orlandi 1994). The present report discusses a set of new 
numerical simulations of isotropic turbulence, and a re-analysis of the old ones, in 
an effort to prove or disprove the presence of these waves in actual turbulent flows 
and to understand the dynamics, as oposed to the kinematics, of the vortices. 

One set of experiments use hyperviscous dissipation instead of regular viscosity. 
Since the strong vortices are known to be dissipative structures with characteris- 
tic radii in the range of the near dissipation range of the turbulence cascade, they 
can be considered as being forced by inertial range eddies, but to be dominated 
by viscous diffusion. The hope was that hyperviscosity would change the diffusion 
mechanism without changing the forcing and give us a better insight into the dy- 
namics. This strategy was largely successful, but some unexpected consequences of 
hyperviscosity were found in the form of a strong bump in the upper inertial range 
of the spectrum, spanning at least one and a half decades in wavenumber, and which 
seems to be an evolution of the much weaker bump described in some experiments 
in real fluids. This bump is discussed below in detail and may be of some practical 
relevance since some sub-grid models in LES computations use eddy viscosity laws 
that resemble hyperviscosity, and the present result suggests that their effect might 
distort the spectrum across a range of scales comparable to that of most practical 
LES numerical grids. 
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Finally, the study of strong vorticity leads directly to the question of intermit- 
tency. It appears that at the Reynolds numbers of our flows (Re \ < 170), most 
of the extended tails in the probability density functions of the velocity gradients 
can be associated to the vortex filaments and scale with Reynolds number in the 
same way as they do. We show evidence that the volume ratio occupied by the fila- 
ments decreases with Re\ while the intensity of the individual structures increases. 
This raises the question, already posed in JWSR, of their behavior in the infinite 
Reynolds number limit. It would appear that the mechanisms described here can 
not be extended indefinitely as Re\ — ► oc although the limitations of numerical 
resolution prevent us from giving a definite answer. 

The simulations, both viscous and hyperviscous, are described first. The global 
flow statistics, especially those related to intermittency, are described next, followed 
by a description of worm properties and a short discussion of the results. The present 
is an interim report, and more work is needed, especially from the point of view of 
theoretical analysis. 

2. Viscous simulations 

The viscous simulations used in this report are essentially the same as in JWSR. 
The spectral numerical scheme (Rogallo 1981) and other parameters are described 
in that paper, which should also be consulted for the detailed definitions of the 
different quantities, which generally follow (Batchelor 1953). Some simulations were 
continued for a longer time to improve the statistics, and the two lower Reynolds 
numbers were repeated at higher resolution, fc max *? = 4 instead of & max 77 = 2, to 
check that the scalings reported in JWSRwere not artifacts of the use of a uniform 
resolution in Kolmogorov units. No artifacts were found. These new simulations 
were also used in a separate study of the effect of resolution on both experiments 
and in simulations (Jimenez 1994b). 

Since the time of the previous report, the fields corresponding to the highest 
Reynolds number (Re \ = 168 on a 512 3 grid) became unavailable for further pro- 
cessing. To paliate this problem, a new simulation was run on a 384 3 grid, resulting 
in Re\ = 141. All these simulations, including the older ones, are sumarised in 
Table 1, which refers only to the numerical resolution k m&x r] = 2. 

As in the original study, all simulations were forced to reach a statistically sta- 
tionary state with a given energy dissipation, adjusted to achieve a desired value 
for k m!tx r]. In particular, every Fourier velocity coefficient with a wave number such 
that k = |k| < 2.5 was multiplied at each time step by a common real factor, chosen 
so that the extra energy introduced in that way was equal to the desired energy 
dissipation. 

It was felt that forcing the flow at such low wave numbers could reduce excessively 
the statistics of the large scales and perhaps affect the flow. To check that effect, 
two simulations were performed in which the same type of forcing was applied to 
the spectral shell 3.5 < |k| <6.5. In general, few differences were found between 
the two types of simulations although the new statistics should clearly be better 
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Re\ 

N 

L< 

LJX 

LJv 

eL/u ' 3 

t/T 

-f 3 

35 

64 

1.8 

2.3 

27 

1.09 

54.2 

0.49 

63 

128 

2.2 

4.2 

65 

0.80 

9.3 

0.50 

94 

256 

2.0 

6.3 

120 

0.72 

8.2 

0.52 

142 

384 

2.4 

9.5 

222 

0.73 

5.9 

0.52 

168 

512 

2.4 

11.2 

286 

0.69 

5.9 

0.52 

26 

128 

0.6 

1.7 

18 

0.96 

11.2 

0.50 

48 

256 

0.7 

3.2 

44 

0.62 

15.1 

0.51 


Table 1. Numerical and flow parameters for the viscous cases analyzed in this 
paper. t/T is the total run time in eddy turnover units, F 3 is the skewness coefficient, 

L is the integral scale, and L t = u' 3 /e is the eddy dissipation scale. The last two 
cases were forced at higher wavenumbers as explained in text. 

than the older ones. Those differences that were found are discussed below in the 
corresponding sections, and the flows themselves are documented in Table 1. The 
better statistics were compensated by the possibility of reaching only much lower 
Reynolds numbers. Fig. 1 presents compensated spectra for the different flows. 

3. Hyperviscosity 

Hyperviscosity, the use of the iterated Laplacian in place of the usual dissipative 
operator, has been used often in the numerical simulation of turbulence in an effort 
to obtain longer inertial ranges out of a given resolution. The hope has been that 
its use would not affect the inviscid aspects of the turbulent flow although it clearly 
changes the characteristics of the dissipative scales. Since the vortices that con- 
cern us seem to be dissipative features, hyperviscous simulations were undertaken 
in the hope of clarifying their dynamics by changing their behavior. Define the 
hyperviscosity exponent a through the momentum equation, 

^ + Vp = (-l)° +1 iA7 2a u (1) 

Dt 

The regular Navier-Stokes equations correspond to a = 1. In spectral calculations, 
the dissipative term is obtained by multiplying the Fourier coefficients of the velocity 
by |k| 2 “, and there is no reason to restrict a to be an integer although the phys- 
ical representation of the dissipation operator may be complicated for non-integer 
exponents. 

The first difficulty in analyzing hyperviscous turbulence simulations is to find 
reasonable scaling quantities. Since the inertial cascade mechanism is not expected 
to change, the Kolmogorov dimensional arguments should apply, and the inner 
microscale should depend only on the viscosity coefficient, v , and of the dissipation, 
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which is expressed in terms of the energy spectrum as 


e = 2u 


/■ 


k la E(k)dk. 


( 2 ) 


The inertial scaling is 

E(k) = v 2 r]f(krj), 

where , , 

v = (ei?) 1 / 3 , 77 = (i/ 3 /e) 1/2(3o-1) . ( 3 ) 

Note that these definitions agree with the usual ones when a = 1, but that v should 
not be interpreted as a regular viscosity coefficient, even dimensionally. 

There is no unique generalization for the Taylor microscale. In regular turbulence 
its most obvious use is in defining the microscale Reynolds number Re x , which is 
then related to the ratio of the different length scales as 

L,ln = 15 , A/ij = (4) 

where L t = u' 3 /e. If we take these two relations as defining A and Re x , the result 
is , /„ 

A = v/l5 (W) 1/3 , i?eA — V^ l ~ a ^ 2a i O) 

which reduces to the usual one for regular viscosity. 

The numerical experiments are listed in Table 2. All except one were conducted 
at the same numerical resolution, 128 3 and k m&x ri = 2, and result in roughly similar 
Re x . A single case was repeated on a larger grid to check Reynolds number scaling. 
All the fields were forced at |k| < 2.5 in the same way as in the viscous cases. One 
of the quantities listed in the table, 

€ * = e/v l/Q u’ (3a - 1)lQ 


is a generalization of the equation e/ uuj' 2 = 1, which only holds for regular viscosity. 
The deviation of this coefficient away from unity measures the failure of the vorticity 
magnitude to represent dissipation in hyperviscous flows. 

S.l Numerics 

There are special numerical problems associated with high hyperviscous expo- 
nents. The principal one is the limitation imposed to the time increment by the 
accuracy requirements of the dissipative term. Since most codes implement this 
substep by some unconditionally stable implicit scheme, the viscous parameter does 
not represent a stability limitation, but unless the time increment is chosen short 
enough, the evolution of the velocity due to dissipation will be represented inaccu- 
rately, and while the resulting system will usually be equivalent to some dissipative 
model, it would not be possible to claim that it represents a hyperviscosity of the 
intended order. 
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a 

•Re a 

N 

L t 

LJ A 

L <h 

eL/u ' 3 

e* 

t/T 

-f 3 

i. 

63 

128 

2.2 

4.2 

65 

0.80 

1.00 

9.3 

0.50 

1.25 

65 

128 

2.2 

4.3 

68 

0.69 

0.57 

3.1 

0.50 

1.50 

64 

128 

2.2 

4.2 

66 

0.68 

0.41 

3.1 

0.49 

2, 

68 

128 

2.4 

4.5 

73 

0.66 

0.27 

3.4 

0.49 

2.50 

68 

128 

2.5 

4.5 

74 

0.65 

0.21 

9.5 

0.49 

4. 

75 

128 

2.8 

5.0 

86 

0.65 

0.15 

4.8 

0.42 

2. 

95 

256 

2.0 

6.3 

121 

0.84 

0.28 

1.5 

0.52 


Table 2. Numerical and flow parameters for the hyperviscous cases analyzed 
in this paper. The case a — 1 is the same one described in Table 1, and is only 
included here for comparison. The quantity e * is defined in the text. 

Consider a spectral code for the evolution equations (1). The accuracy and sta- 
bility of the convective term is controlled by the usual CFL parameter, which is 
proportional to 

c ~ u fc max A*, (6) 

while that of the hyperviscous term is controlled by 

6 = (7) 

In our code, which uses an integrating factor for the dissipative term (Rogallo 1981), 
the modes corresponding to fc max are multiplied every time step by exp (—(5), and 
the explicit integrator for the advection term assumes that this factor is close enough 
to unity for the accuracy of the nonlinear interactions not to be destroyed during 
the step. This requires that 6 should not be too large. Other implicit schemes have 
different behaviors, but the accuracy requirement for the dissipative term is always 
that 6 be at most 0(1). 

The time increment At is usually adjusted by fixing the CFL, and the viscous 
parameter has to be measured. In our code, u 1 is given by the dissipation ( and 
by the form of the spectrum through the choice of the resolution parameter k mAX ri. 
Assume that e ~ u' 3 / L ~ u’ 3 , where l'o is the center of our forcing band. It 
follows from our definition of the Kolmogorov scale that 

u ' ~ iv -2 "(fc 0 f?r i/3 . (8) 

When At is obtained using this estimate in the CFL and is then substituted in 6, 
we obtain 

6/c = const, (fcma xf?) 20 ~ 2/3 (fcoA mi «) 1/3 « 1- 
From our observations, the empirical constant is approximately 0.03. 


( 9 ) 
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Figure 2. a) Dissipation spectra, e _2 / 3 r? _5/3 (fc?j) 2Qr ^(A:), for different hypervis- 
cosity exponents. From left to right: a = 1, 1.25, 1.5, 2.0, 2.5, 4.0. b) Wavenumber 
of maximum dissipation. 

Fig. 2a displays the hyperviscous dissipation spectra k 2a E(k ), for various expo- 
nents, and Fig. 2b shows the position of the dissipation maximum for the different 
cases. It is clear from these figures that hyperviscosity succeeds in separating the 
forcing and dissipation ranges, even at moderate resolutions, although it will be 
seen below that this does not guarantee an inertial range in the usual sense. It is 
also clear that fc max 7 / should be chosen somewhat larger than unity if the dissipa- 
tive range is to be reasonably well represented, making the requirement expressed 
in Eq. (9) more restrictive as the hyperviscosity order increases. In practice, for our 
choice of fc max r ? = 2 and for fc 0 - 2, fc max = 60 - 120 (128 3 or 256 3 grids), it limits 
our simulations to a < 4 — 5. 

3.2 The near dissipation spectrum 

Compensated three-dimensional spectra for some of the hyperviscous fields are 
presented in Fig. 3a. There is no collapse in the dissipative range, and the quali- 
tative shape of the spectra changes drastically as a function of the hyperviscosity 
exponent although it should be noted that the two spectra corresponding to the 
same exponent, a = 2, collapse well. As the exponent increases, a large “bump” 
appears in the near-dissipation range that for the highest exponents and at these 
low Reynolds numbers, dominates the spectrum. In Fig. 3b, in which the spectrum 
is premultiplied by fc, the bump is seen to behave approximately as k~ x . 

Anomalous bumps in this region of the spectrum have been reported in experi- 
mental viscous flows (Mestayer 1982, She & Jackson 1993, Saddoughi &; Veeravalli 
1994), and the fc” 1 behavior was claimed by the second of these groups. A theo- 
retical explanation for this power law, based on the depletion of nonlinearity in the 
near dissipation scales, was offered in (Yakhot & Sakharov 1993). An older expla- 
nation for the presence of the bump, although not for the power behavior, is that 
the outgoing energy cascade is inhibited in the dissipation range since energy has to 
move into spectral triads with much smaller amplitudes than those corresponding 
to the equilibrium k~ 5//3 spectrum while the incoming energy is not subject to that 
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FIGURE 4. Compensated transverse spectra showing extent of near dissipation 

bump. : Rc\ = 168, ol = 1; dashed: Re\ ~ 64, a = 1.5; : Re\ = 

96, a = 2; : Rt\ = 64, a = 4; : Re\ = 600, a = 1, after (Saddoughi k 

Veeravalli 1994). The peak in the lowest wave numbers of all the numerical spectra 
is spurious due to forcing scheme. 


inhibition. The bump is a consequence of the resulting enery “bottleneck” in the 
cross-over region (Falkovich 1994). The latter author has described the same spec- 
tral behavior in wave turbulence (Ryzhenkova k Falkovich 1990) and has remarked 
that it would become more pronounced in hyperviscous simulations as the damping 
becomes more abrupt. The effect is present in spectral turbulence models and has 
been observed in EDQNM based simulations (Mestayer, Chollet k Lesieur 1983). 
Finally, some recent hyperviscous simulations of isotropic turbulence, using a = 8, 
have also reported the presence of a A:*” 1 bump (Borue k Orszag 1994). 

Since it is known that the near dissipation region contains strong vortex filaments, 
they may also provide an explanation for the fc -1 behavior. That randomly oriented 
vortex filaments would generate a spectrum with this behavior was first noted by 
Townsend (1951) and observed directly in JWSR by computing the spectrum of 
a flow in which all the vorticity, except the one contained in the strong filaments, 
had been zeroed. Moreover, an (inviscid) infinite vortex filament is an equilibrium 
solution of the Euler equation for which, by definition, all the nonlinear interactions 
cancel identically, and for which the turbulent cascade is fully absent, even if the 
filament itself is usually formed by the cascade mechanism of vortex stretching. In 
this sense the three explanations are not necessarily incompatible, and the filament 
hypothesis merely points to a possible physical implementation of the two previous 
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ideas. 

However, it will be shown below that in the course of the present investigation, 
we were not able to find significant differences between the structure or frequency 
of filaments in viscous and hyper viscous flows, which may explain the presence of 
stronger bumps in the latter. As a consequence the filament explanation remains 
unlikely. Also, our observations cast doubts on the —1 exponent as a preferred 
spectral law. 

In Fig. 4 we present transverse compensated spectra from several simulations, 
including our highest Reynolds number regular viscosity case, together with a ex- 
perimental spectrum from (Saddoughi Sz Veeravalli 1994). The bump is clearly seen 
in the experimental spectrum and, in retrospect, also in the regular viscous simu- 
lation, showing that even the highest Reynolds number computed by us is far from 
containing a true inertial range. 

Figs. 5a-b display the local logarithmic slope, d(log E)/d(log k), for the different 
spectra in the near dissipation region. They also include data from (Saddoughi 
& Veeravalli 1994) although in this case, to make it comparable to ours, it is in 
the form of a three-dimensional spectrum obtained from their data for En, using 
the assumption of isotropy. Fig. 5a, for the regular viscous cases, reinforces the 
conclusion that all our simulations are far from the real inertial range. Fig. 5 b, which 
contains the hyperviscous cases, shows a steady climb of the least negative slope 
with increasing a with no sign of saturating at Moreover, the two Reynolds 

numbers at a = 2 show the same tendency as the regular viscous cases with the 
spectral hump extending further into the low wavenumber range, suggesting that 
these simulations are also far from reaching a true inertial range and that the extent 
of the bottleneck region is not shortened by the sharpening of the dissipation peak. 

4. Velocity gradients 

It has been known for a long time that the probability density functions for 
velocity differences become increasingly non-Gaussian at short distances, and that 
this leads to highly intermittent behavior of the velocity gradients, which becomes 
more pronounced at high Reynolds numbers. In JWSRwe presented histograms 
for different combinations of gradients (vorticity and total strain magnitudes and 
local stretching, a = u>Su>/u; 2 ), which clearly showed this Reynolds number effect. 
A compact representation of this departure from Gaussianity is provided by the 
high order flatness and skewness -Fn(C) =< ( n > / < ( 2 > n ^ 2 , two of which axe 
given in Fig. 6 for the viscous and hyperviscous cases. In each case the vertical 
origin of coordinates has been chosen to coincide with the Gaussian value (F 4 = 3 , 
F$ = 15. The tendency for both flatnesses to grow with Reynolds number is clear, 
as is the fact that the simulations forced at intermediate wave numbers do not 
differ significantly in this respect from those forced at low wave numbers. The 
values for the fourth order flatness of the longitudinal gradients agrees well with 
the compilation in (Van Atta &: Antonia 1980) although our experiments cover a 
much narrower range of Reynolds numbers than theirs. 

The variation of the flatness with the hyperviscosity exponent (Fig. 6 b) is more 
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Figure 5. Logarithmic slope for the different three-dimensional energy spectra, 

in the near dissipation region, a) Regular viscosity. : Re\ = 168; : 141; 

: 98; : 65. b) Hyperviscosity. : a = 1.25; : 1.5; : 

2. (Re\ ~ 68); +: 2. (Re \ = 94); : a = 4. In both figures, open circles are 

experimental data from (Saddoughi &: Veeravalli 1994) at Re\ = 600. 
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A A A 



FIGURE 6 . Fourth and sixth order flatness factors for velocity gradients, (a) Vis- 
cous simulations vs. Reynolds number, (b) Hyperviscous simulations vs. exponent; 
Re\ « 65. ^ : longitudinal gradients, du/dx. o : transverse, du/dy. Open symbols 
are flows forced at low wave numbers. Closed symbols are forced at intermediate 
wave numbers. 


surprising. While it seems from table 2 that the equivalent Reynolds number in- 
creases slightly with increasing a and while the same conclusion could be reached 
from the increasing separation between the integral scales and the dissipation peak 
(Fig. 2a), the flatness is seen to decrease with a. The same behavior is observed in 
the histograms of almost all the gradient quantities, as can be seen for the vorticity 
magnitude in Fig. 7a. The exception is the stretching histogram (Fig. 7 b) which 
shows a much weaker variation with a, essentially within the statistical uncertainty. 

It should be stressed that the variation of the flatness with .Re* is similar in 
the hyperviscous flows at constant a and in the viscous ones. This is clear from 
table 3, which contains flatness for the two Re\ available at a = 2 . There is also 
a consistent, although much weaker, increase in the skewness, F 3 , with Re \ which 
can be seen from tables 1 and 2 , and which is also in general agreement with the 
data in (Van Atta & Antonia 1980). 

5. Worms 

An algorithm was developed in JWSR to track individual filaments and to mea- 
sure their properties. The same algorithm has been used here to analyze the new 
flow fields and to obtain data which may be useful in understanding their dynamics. 
The axis of each filament is followed until some arbitrary definition of its end point 
is reached and the vortex radius, circulation, peak axial vorticity, and axial stretch- 
ing are measured at each point and used to compile statistics over the fraction of 
the flow contained in the worms. 

The axial stretching is defined as , where Sij is the strain tensor and rii 

is the unit vector in the direction of the axis. The condition to end a filament was 
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Re x 

F±i 

f 4 , 

Foe 

Fet 

68 

4.5 

5.9 

49 

100 

95 

4.8 

6.5 

62 

131 


TABLE 3. Flatness factor for the two hyperviscous cases with a = 2. Subscript £ 
refers to longitudinal gradients, t to transverse ones. 

described in detail in JWSR but depends essentially on the axial vorticity weakening 
below a given level. 

The simplest statistic is the mean value of a given property over all the axial 
points of all the worms detected in a given flow. There are two groups of properties. 
Since filaments have to be formed by stretching of preexisting vorticity, the axial 
stretching can be considered as the driving force in the formation of the worms 
while all the other properties can be considered as the results of that forcing. 

Both groups behave differently. The mean value of axial stretching, averaged 
over the axes of the worms, is shown in Fig. 8 normalized with the rms value 
of the vorticity over the whole flow, which is the natural scaling for the velocity 
gradients in the bulk of the flow. Except for a slight tendency to grow with Re\, 
the stretching scales well with u/ although the proportionality constant depends on 
the hyperviscosity exponent. It was shown in JWSR that the average stretching, 
conditioned on a given vorticity, depends only weakly on the conditioning value, 
and this was interpreted as an indication that the filament did not stretch itself 
appreciably. The present result reinforced that conclusion since it will be seen 
below that the mean vorticity in the worms is a fairly strong function of Re\ while 
we have shown that the stretching is not. 

This is further confirmed by Figs. 9 and 10, which show correlation functions and 
lengths for the axial distribution of the stretching along worms. The correlation 
length is a useful indication of the typical length scale of a given property, and we 
use it here as a less subjective indication of length than the one used in JWSR, 
which was based on the total length of the detected worms. 

Consider a quantity, £(s), given as a function of axe length along the axis and 
assumed to have zero mean value, and define the correlation 

R(s) =< C(a')C(a' + s) > / < C 2 >, (10) 

where the average is taken over the axial position s' . This function has zero integral 
and i?(0) = 1. It will first become negative at some point s = sq (Fig. 9). We define 
the correlation length as 

= f sR(s)ds/ f R(s)ds. (11) 

J 0 Jo 

The same definition can be used on properties which are defined at all points in 
the flow for which the arc length s has to be substituted by a coordinate distance 
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FIGURE 8. Average axial stretching along filament axes. : viscous flows; o : 

a = 2; * : a = 4. 



s/t) 


FIGURE 9. Correlation functions of axial stretching, against of the arc length 

along the worm axes. : a = 1 , Re\ = 63; : Rz\ = 142; : 

a = 2, Re\ — 68; : a = 4. 
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Figure 10. Axial correlation length of the stretching along the filament axes, 
(a) Normalized with the Kolmogorov length, (b) Normalized with the correlation 
length of longitudinal gradients over the whole flow. Symbols as in Fig. 8. The 
dashed line in (a) is the Taylor microscale, for comparison. 

(e.g. x). This global correlation function can be computed as the inverse Fourier 
transform of the one-dimensional spectrum of the quantity (Batchelor 1953). From 
that, and since the power spectrum of the velocity gradients has its maximum at 
the near disipation region, it follows that the correlation lengths for gradients when 
averaged over the whole flow should scale with 77. 

Fig. 10a shows the correlation length for the stretching along the worms, £„, which 
scales well with r), although the scaling constant depends again on the hyperviscosity 
exponent. This is not surprising since we have shown that the shape of the near 
dissipation spectra depends on a and so does presumably the relation between the 
size of the gradient eddies and 77. In Fig. 10b we show ( a normalized with the 
correlation length of the longitudinal gradients du/dx taken over the whole flow. 
Most of the dependence with a disappears, and the ratio is close to one, showing 
once again that the stretching along the worms is essentially indistinguishable from 
that at a generic point in the flow and that self stretching is not important. 

It is therefore surprising that the axial correlation lengths of all the resulting 
worm properties are much larger than 7/ and scale apparently on the integral scale 
of the flow (Fig. 11). This generalizes the observation in JWSR, which has been 
made by all the investigators that have studied these filaments, that their length 
is of the order of the integral scale. While the original observations were made on 
the basis of arbitrary visualization or termination criteria, the present result is a 
more objective characterization of the same phenomenon. In these last two figures 
we have included a curve for the Taylor microscale A which is defined by Eq. (4). 
This is done mainly for comparison, but also because this scale has been suggested 
often as a natural scale for the turbulent eddies. In particular, it was observed 
in JWSR that if the characteristic stretching along worms was of 0(w ') and since 
it is known that the velocity differences in a turbulent flow scale with u' without 
strong Reynolds number effects, the largest possible length scale for the stretching 
was 0(u>' /u') = 0(A). It was felt at the time that this contradicted the observation 
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FIGURE 11. Axial correlation lengths along worms axes, normalized of the dissi- 
pation length L ( . and o : axial vorticity at axis; and x: circulation; 

and o : radius. Lines without symbols are viscous simulations. Lines with 

symbols axe axe a = 2. Isolated symbols are a = 4. 



FIGURE 12. Probability density function for the stretching a along the axes of the 
worms. Symbols as in Fig. 9. All histograms are normalized to unit mean value. 
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Figure 13. Mean worm radius, (a) Normalized with the Kolmogorov length, (b) 
Normalized with Burgers’ radius for the average axial strain. Symbols as in Fig. 8. 


that the length of the worms was 0(L ( ) since it was assumed that vortices could 
only survive if stretched on the average. This proved to be false, and it was shown 
in (Jimenez & Wray 1994) that the strain along the worms is compressive over a 
large fraction of the axis (Fig. 12). Verzicco, Jimenez and Orlandi (1994) have 
shown recently by direct computation that vortices subject to spatially nonuniform 
axial strain, even one with zero axial average stretching, can reach a steady state 
thanks to the presence of axial inertial waves which “smooth” the strain over the 
compressive parts. 

It is still interesting to note that although the L t scaling seems to be a better 
representation of the data in Fig. 11, especially at the higher Re\, a scaling of a 
few Taylor microscales is not completely incompatible with them. 

In (Jimenez &: Wray 1994) the equivalent to Fig. 12 was plotted together with 
probability density functions of the same quantity computed over the whole flow 
field. The differences between the two were shown to be small, stressing again that 
the forcing of the filaments is not different from that of the rest of the flow. Note 
however that the collapse of the stretching pdf’s over the full field for different values 
of a (Fig. 7b) is better than those of the same quantity over the worms (Fig. 12), 
suggesting either that there is a selection mechanism for the location of worms in 
terms of stretching or a dynamical feedback from the worms into the structure of 
the stretching itself. The same conclusion can be drawn from the dependence of the 
average worm stretching on a, shown in Fig. 10. 

The average radius of the worms is given in Fig. 13a. It scales with the Kol- 
mogorov scale although there is again a different proportionality constant for dif- 
ferent hyperviscosity exponents. A better collapse is possible if we assume that the 
filaments are Burgers’ vortices driven by the mean axial strain. For viscous flows, 
the Burgers’ radius due to a strain u/ is 2(u/u>') 1 / 2 = 2 rj, which agrees with the 
approximate Kolmogorov scaling for these flows in Fig. 13a. 

Hyperviscous Burgers’ vortices were computed in (Jimenez 1994a). Although 
their outer tails are different from the viscous ones and actually change sign before 





Small-scale vorticity in isotropic turbulence 


305 




FIGURE 14. (a) Radial vorticity distribution for hyperviscous Burgers’ vortices, 

normalized using Eq. (12). : a — 1; : 2; : 3; : 4. (adapted 

from Jimenez 1994a). (b) Mean radial distribution of axial vorticity in worms. 

Symbols as in Fig. 9. 
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they settle to zero, their cores are approximately Gaussian (Fig. 14a) with a 1/e 
radius which can be approximated for a < 5 by 

R b = (2.022 - 0.0843a + 0.0688a 2 - 0.0086a 3 )(i//a) 1 / 2a , (12) 

where the coefficient is the result of a polynomial fit to the solution of an ordinary 
differential equation. When this formula is used to compute the Burgers’ radius as- 
sociated to the measured mean stretching from Fig. 8, the results are the normalized 
radii in Fig. 13b, in which the effect of hyperviscosity is essentially absent. 

Fig. 14. b contains vorticity profiles averaged along the axes of all the worms in 
a given flow. The profiles are also approximately Gaussian, and it is interesting 
to note that the more hyperviscous profiles show signs of negative vorticity at the 
edges, which may be a reflection of the oscillations in the tails of the theoretical 
solutions. 

Full probability density functions for the local radius and for the radius divided 
by the local Burgers’ radius are given in Fig. 15. To avoid effects related to the 
variations of the mean values, the averages of all those distributions are adjusted to 
unity. The collapse of the Burgers’ plot is excellent, but there is some differences 
between the pdf’s of the raw radii in viscous and in hyperviscous flows, probably 
reflecting the differences in the structure of the local stretching that were discussed 
above. 

The mean values of the vorticity at the axes of the worms are given in Fig. 16a. 
It was suggested in JWSR that u>/u/ increases with Reynolds number as Re]/ 2 and 
that same normalization is used in the figure. Except for the case at Re\ = 48 
in which the forcing was done at higher wave numbers, it seems to work correctly. 
The intermediate forcing seems to work differently from all the other flows for all 
quantities which scale with the integral length as can also be observed in Fig. 11. 
The large scales are different from those of the flows forced at lower wave numbers, 
and it appears that the axial distribution of vorticity in the worms is controlled 
by them. A line corresponding to uj ~ u/ is also included in the figure and is 
incompatible with the data. In Fig. 16b we have represented the mean value of the 
vortex circulation, normalized with the mean vorticity and radius for each flow. It 
clusters around unity in what is essentially a consistency check for the averages, 
but which also reflects that the statistical distribution of vorticity and radius are 
relatively independent of the dissipation model. A few representative pdf’s are 
given in Fig. 17, normalized to unit mean. The sharp cut-off of the vorticity pdf 
is artificial. The tracking algorithm terminates a worm whenever its axial vorticity 
falls below u/ . 

It is finally interesting to enquire which is the relative importance of worms 
with respect to the bulk of the flow. This is largely a matter of definition, but a 
volume fraction can be defined by taking the mean vorticity detected at the axes 
by our algorithm and defining as worms all the points whose vorticity magnitude is 
above that threshold. This, although arbitrary, seems justified since it appears from 
Fig. 17a that our threshold is below the maximum of the distribution and is probably 
not distorting the mean value too much. When this is done, a volume fraction can 
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Figure 16. (a) Mean vorticity along worm axes, normalized with u'Re] /2 . Sym- 

bols as in Fig. 8 . Dashed line corresponds to u> ~ w', without Reynolds number 
dependence, (b) Mean value of filament circulation, normalized with mean vorticity 
and vortex radius, < 7 > /n < u; >< r > 2 . 


be read from the vorticity histogram of the flow (Fig. 7a) and is represented in 
Fig. 18. The volume occupied by the worms decreases with Reynolds number as 
Re \ 2 - Since their characteristic vorticity increases only with Re]/ 2 , this implies that 
not only the volume, but also the enstrophy contained in them decreases rapidly as 
the Reynolds number increases. 

6. Conclusions and future work 

We have shown that the vortex strong filaments in isotropic turbulence have 
lengths and axial correlation lengths of the order of the integral scale of the flow, 
even if the stretching along their axes has a spatial scale of the order of the Kol- 
mogorov length and seems essentially indistinguishable from the strain in the bulk 
of the flow. We have also shown that the average radius of the vortices is very 
close to the Burgers’ radius corresponding to the mean axial strain, even if large 
segments of the axes are actually under compression. We have suggested that this is 
accomplished through the action of axial waves, which distribute the effective strain 
along the axes. All these observations hold for hyperviscous flows with the obvious 
modifications needed to accommodate the different core structure of the vortices. 

We can think of few mechanisms to generate coherent vortices of such lengths 
in turbulent flows. The obvious one, in which vortices form by roll-up of the large 
scale vorticity layers in between large eddies, is unlikely because it can readily be 
shown that the large-scale stretching is at most able to collapse vortices to radii 
of the order of the Taylor microscale, and it is difficult to think of a way in which 
straining motions of the observed scales, r), could further collapse these cores, much 
thicker than themselves, into more compact cores. 

Another possibility, and the one that we favor at the moment, is a mechanism 
by which short vortex “sticks” form individually and are later patched by the axial 
waves into longer units. This mechanism has been demonstrated in simpler situa- 
tions in (Verzicco, Jimenez and Orlandi 1994). Note that from Figs. 10 and 11 , the 
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FIGURE 17. Probability density functions for (a) vorticity at worm axis, (b) 
Circulation, adjusted in each case to unit average. Symbols as in Fig. 9. 
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Figure 18. Volume fraction occupied by vorticity above the average detected in 
worms, o : viscous simulations, low wave number forcing; • : viscous, middle wave 
number forcing; a : hyperviscous; dashed line has slope -2. 


elongation of the worms, i r / < r >, although large and increasing with Re\, never 
exceeds 10-20 in our range of Reynolds numbers. 


A remaining problem is the growth of the circulation of the vortices with Re\. 
It was noted in JWSR that this implies an increase of the Reynolds number of 
the azymuthal motion of the individual structures and that at some point those 
structures should become unstable. We still feel that way, but we have no way of 
predicting the conditions for this instability. 


We have also shown that the volume fraction occupied by worms decreases as 
This is somewhat misleading and should be understood in relation to the 


Re? 


decrease of their radii. If we remember that rj ~ Re x 3/2 , it follows that the accu- 
mulated length of the structures grows like L ( Re\, and their total “number” per 
unit volume, as Re\. This is difficult to reconcile with their increasing instability 
as their Reynolds number increases. It has been proposed recently that an internal 
structural transition might exist in turbulence at Re\ ss 10 3 , resulting in a decrease 
of intermittency (Tabeling et al. 1994). If this is confirmed, it could possibly resolve 
the present dilemma. 


Finally, we have described the effect of hyperviscosity in isotropic simulations. 
The resulting spectra are dominated by large humps in the near dissipation region, 
which extend to the last wave number decade of the inertial range. They seem to 
be generated by the energy “bottleneck” produced by the inhibition of the energy 
cascade by viscosity. 
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